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Cosmic rays are particles (mostly protons) accelerated to relativistic speeds. 
Despite wide agreement that supernova remnants (SNRs) are the sources of 
galactic cosmic rays, unequivocal evidence for the acceleration of protons in 
these objects is still lacking. When accelerated protons encounter interstellar 
material they produce neutral pions, which in turn decay into gamma rays. 
This offers a compelling way to detect the acceleration sites of protons. The 
identification of pion-decay gamma rays has been difficult because high-energy 
electrons also produce gamma rays via bremsstrahlung and inverse Compton 
scattering. We detected the characteristic pion-decay feature in the gamma- 
ray spectra of two SNRs, IC 443 and W44, with the Fermi Large Area Tele- 
scope. This detection provides direct evidence that cosmic-ray protons are 
accelerated in SNRs. 

A supernova explosion drives its progenitor material supersonically into interstellar space, 
forming a coUisionless shock wave ahead of the stellar ejecta. The huge amount of kinetic en- 
ergy released by a supernova, typically 10^^ ergs, is initially carried by the expanding ejecta 
and is then transferred to kinetic and thermal energies of shocked interstellar gas and relativistic 
particles. The shocked gas and relativistic particles produce the thermal and nonthermal emis- 
sions of a supernova remnant (SNR). The mechanism of diffusive shock acceleration (DSA) 
can explain the production of relativistic particles in SNRs (i). DSA generally predicts that a 
substantial fraction of the shock energy is transferred to relativistic protons. Indeed, if SNRs 
are the main sites of acceleration of the galactic cosmic rays, then 3 to 30% of the supernova ki- 
netic energy must end up transferred to relativistic protons. However, the presence of relativistic 
protons in SNRs has been mostly inferred from indirect arguments (2-5). 

A direct signature of high energy protons is provided by gamma rays generated in the 
decay of neutral pions (7r°); proton-proton (more generally nuclear-nuclear) collisions create 
7r° mesons which usually quickly decay into two gamma rays (6-8) (schematically written as 
p + p — 7> other products, followed by 7r° — > 27), each having an energy of mjr()C^/2 = 67.5 
MeV in the rest frame of the neutral pion (where rriT^o is the rest mass of the neutral pion and c is 
the speed of light). The gamma-ray number spectrum, F{e), is thus symmetric about 67.5 MeV 
in a log-log representation (9). The 7r°-decay spectrum in the usual e'^F{e) representation rises 
steeply below ~ 200 MeV and approximately traces the energy distribution of parent protons 
at energies greater than a few GeV This characteristic spectral feature (often referred to as the 
"pion-decay bump") uniquely identifies 7r°-decay gamma rays and thereby high-energy protons, 
allowing a measurement of the source spectrum of cosmic rays. 

Massive stars are short-lived and end their lives with core-collapse supernova explosions. 
These explosions typically occur in the vicinity of molecular clouds with which they interact. 
When cosmic-ray protons accelerated by SNRs penetrate into high density clouds, 7r°-decay 
gamma-ray emission is expected to be enhanced because of more frequent pp interactions rel- 
ative to the interstellar medium (10). Indeed, SNRs interacting with molecular clouds are the 
most luminous SNRs in gamma rays (11, 12). The best examples of SNR-cloud interactions in 



our galaxy are the SNRs IC 443 and W44 (13), which are the two highest-significance SNRs 
in the second Fermi Large Area Telescope (LAT) catalog (2FGL) {14) and are thus particularly 
suited for a dedicated study of the details of their gamma-ray spectra. The age of each remnant 
is estimated to be ~ 10, 000 years. IC 443 and W44 are located at distances of 1.5 kpc and 2.9 
kpc, respectively. 

We report here on 4 years of observations with the Fermi LAT (4 August 2008 to 16 July 
2012) of IC 443 and W44, focusing on the sub-GeV part of the gamma-ray spectrum - a cru- 
cial spectral window for distinguishing vr^-decay gamma rays from electron bremsstrahlung or 
inverse Compton scattering produced by relativistic electrons. Previous analyses of IC 443 and 
W44 used only 1 year of Fermi LAT data (15-17) and were limited to the energy band above 
200 MeV, mainly because of the small and rapidly changing LAT effective area at low energies. 
A recent update to the event classification and background rejection (so-called Pass 7) provides 
an increase in LAT effective area at 100 MeV by a factor of ~ 5 (18), enabling the study of 
bright, steady sources in the galactic plane below 200 MeV with the Ferm?'-LAT. Note that the 
gamma-ray spectral energy distribution of W44 measured recently by the AGILE satellite falls 
steeply below 1 GeV, which the authors interpreted as a clear indication for the 7r°-decay origin 
of the gamma-ray emission (19). Also, a recent analysis of W44 at high energies (above 2 GeV) 
has been reported (20), revealing large-scale gamma-ray emission attributable to high-energy 
protons that have escaped from W44. Here we present analyses of the gamma-ray emission 
from the compact regions delineated by the radio continuum emission of IC 443 and W44. 

The analysis was performed using the Fermi LAT Science Tools {21). We used a maxi- 
mum likelihood technique to determine the significance of a source over the background and 
to fit spectral parameters (22, 23). For both SNRs, additional sources seen as excesses in the 
background-subtracted map have been added to the background model (24) and are shown as 
diamonds in Fig. [T]- one in the case of IC 443, three in the case of W44. The inclusion of 
these sources had no influence on the fitted spectrum of the SNRs. Three close-by sources 
(2FGL J1852.8-I-0156C, 2FGL J1857.2-i-0055c, and 2FGL J1858.5-i-0129c) have been identified 
with escaping cosmic rays from W44 {20). These 2FGL sources have been removed from the 
background model (see below) in order to measure the full cosmic-ray content of W44. 

Figure |2] shows the spectral energy distribution obtained for IC 443 and W44 through max- 
imum likelihood estimation. To derive the flux points we performed a maximum likelihood 
fitting in 24 independent logarithmically spaced energy bands from 60 MeV to 100 GeV. The 
normalization of the fluxes of IC 443 and W44 and those of neighboring sources and of the 
galactic diffuse model, was left free in the fit for each bin. In both sources, the spectra below 
~ 200 MeV are steeply rising, clearly exhibiting a break at ~ 200 MeV. To quantify the sig- 
nificances of the spectral breaks, we fit the fluxes of IC 443 and W44 between 60 MeV and 2 
GeV - below the high-energy breaks previously found in the 1-year spectra {15, 16) with both 
a single power law of the form F{e) = K{e/eo)~^^ and a smoothly broken power law of the 
form F{e) = K{e/eo)~^'il + (e/^br)^^'"^'^/")"" with Eq = 200 MeV. The spectral index 
changes from Fi to r2 (> Fi) at the break energy ^br- The smoothness of the break is deter- 
mined by the parameter a, which was fixed at 0.1 (Table [I])- We define the test-statistic value 



(TS) as 21n(£i/£o) where £1/0 corresponds to the likelihood value for the source/no- source 
hypothesis (23). The detection significance is given by ~ y/TS. The smoothly broken power 
law model yields a significantly larger TS than a single power law, establishing the existence of 
a low-energy break. The improvement in log likelihood when comparing the broken power law 
to a single power law corresponds to a formal statistical significance of 19cr for the low-energy 
break in IC 443 and 21a for that in W44, when assuming a nested model with two additional 
degrees of freedom. 

Table 1 : Spectral parameters in the energy range of 60 MeV to 2 GeV for power-law (PL) and 
broken power-law (BPL) models. TS = 2ln{Ci/ C2) is the test-statistic value. 
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To determine whether the spectral shape could indeed be modeled with accelerated protons, 
we fit the LAT spectral points with a 7r°-decay spectral model, which was numerically calculated 
from a parameterized energy distribution of relativistic protons. Following previous studies 
{15, 16), the parent proton spectrum as a function of momentum p was parameterized by a 
smoothly broken power law in the form of 



dK 



dp 



p_ 

Phi 



(1) 



Best-fit parameters were searched using x^-fitting to the flux points. The measured gamma-ray 
spectra, in particular the low-energy parts, matched the 7r°-decay model (Figure |2l). Parameters 
for the underlying proton spectrum are si = 2.36 ± 0.02, S2 = 3.1 ± 0.1, and = 239 ± 
74 GeVc~i for IC 443 and si = 2.36 ± 0.05, Sa = 3.5 ± 0.3, and Pbr = 22 ± 8 GeY c'^ for 
W44 (statistical errors only). In Figure |3] we show the energy distributions of the high-energy 
protons derived from the gamma-ray fits. The break pbr is at higher energies and is unrelated to 
the low-energy pion-decay bump seen in the gamma-ray spectrum. If the interaction between 
a cosmic-ray precursor (i.e., cosmic rays distributed in the shock upstream on scales smaller 
than ~ 0.1-R, where R is the SNR radius) and adjacent molecular clouds were responsible for 
the bulk of the observed GeV gamma rays, one would expect a much harder energy spectrum 
at low energies (i.e. a smaller value for the index si), contrary to the Fermi observations. 
Presumably, cosmic raysin the shock downstream produce the observed gamma rays; the first 



index si represents the shock- acceleration index with possible effects due to energy-dependent 
propagation, and pbr may indicate the momentum above which protons cannot be effectively 
confined within the SNR shell. Note that pbr results in the high-energy break in the gamma-ray 
spectra at ~ 20 GeV and ~ 2 GeV for IC 443 and W44, respectively. 

The 7r°-decay gamma rays are likely emitted through interactions between "crushed cloud" 
gas and relativistic protons, both of which are highly compressed by radiative shocks driven into 
molecular clouds that are overtaken by the blast wave of the SNR (25). Filamentary structures 
of synchrotron radiation seen in a high-resolution radio continuum map of W44 (26) support 
this picture. High-energy particles in the "crushed cloud" can be explained by re-acceleration of 
the pre-existing galactic cosmic rays (25) and/or freshly accelerated particles that have entered 
the dense region (20). The mass of the shocked gas (~ 1 x IO^Mq and ~ 5 x IO^Mq for IC 443 
and W44 respectively, where Mq is the mass of the Sun) is large enough to explain the observed 
gamma-ray luminosity. Because the "crushed cloud" is geometrically thin, multi-GeV particles 
are prone to escape from the dense gas, which may explain the break pbr- 

Escaped CRs reaching the unshocked molecular clouds ahead of the SNR shock can also 
produce 7r°-decay gamma rays (27, 28). Indeed, the gamma rays emitted by the escaped CRs in 
the large molecular complex that surrounds W44 (total extent of 100 pc) have been identified 
with three close-by sources (20), which led us to remove them from the model in the maximum 
likelihood analysis, as mentioned above. With this treatment, the measured fluxes below 1 GeV 
contain small contributions from the escaped CRs, but this does not affect our conclusions. The 
escaped CRs may significantly contribute to the measured TeV fluxes from IC 443 (29, 30). 
Emission models could be more complicated. For example, the CR precursor with a scale 
of ~ 0.1-R at the highest energy could interact with the adjacent unshocked molecular gas, 
producing hard gamma-ray emission. This effect is expected to become important above the 
LAT energy range. 

We should emphasize that radiation by relativistic electrons can not as naturally account for 
the gamma-ray spectra of the SNRs (24). An inverse-Compton origin of the emission was not 
plausible on energetic grounds (11). The most important seed photon population for scattering is 
the infrared radiation produced locally by the SNR itself with an energy density of ~ 1 eV cm~^, 
but this is not large enough to explain the observed gamma-ray emission. Unless we introduce 
in an ad hoc way an additional abrupt break in the electron spectrum at 300 MeV (Fig. [2] 
dash-dotted lines), the bremsstrahlung models do not fit the observed gamma-ray spectra. If 
we assume that the same electrons are responsible for the observed synchrotron radiation in the 
radio band, a low-energy break is not expected to be very strong in the radio spectrum and thus 
the existing data do not rule out this scenario. The introduction of the low-energy break intro- 
duces additional complexity and therefore a bremsstrahlung origin is not preferred. Although 
most of the gamma-ray emission from these SNRs is due to 7r°-decay, electron bremsstrahlung 
may still contribute at a lower level. The Fermi LAT data allow an electron-to-proton ratio Kep 
of ~ 0.01 or less, where Kep is defined as the ratio of dNe/dp and dNp/dp at p = 1 GeVc~^ 
(figs.|S2]and[S3]). 

Finding evidence for the acceleration of protons has long been a key issue in attempts to elu- 



cidate the origin of cosmic rays. Our spectral measurements down to 60 MeV enable the iden- 
tification of the 7r°-decay feature, thus providing direct evidence for the acceleration of protons 
in SNRs. The proton momentum distributions, well-constrained by the observed gamma-ray 
spectra, are yet to be understood in terms of acceleration and escape processes of high-energy 
particles. 
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Figure 1 : Gamma-ray count maps of the 20° x 20° fields around IC 443 (left panel) and W44 
(right panel) in the energy range 60 MeV to 2 GeV. Nearby gamma-ray sources are marked as 
crosses and squares. Diamonds denote previously undetected sources. For sources indicated 
by crosses and diamonds, the fluxes were left as free parameters in the analysis. Events were 
spatially binned in regions of side length 0.1°, the units of the color bar is the square root of 
count density, and the colors have been clipped at 20 counts per pixel to make the galactic 
diffuse emission less prominent. Given the spectra of the sources and the effective area of the 
LAT instrument, the bulk of the photons seen in this plot have energies between 300 and 500 
MeV. IC 443 is located in the galactic anti-center region, where the background gamma-ray 
emission produced by the pool of galactic cosmic rays interacting with interstellar gas is rather 
weak relative to the region around W44. The two dominant sources in the IC 443 field are the 
Geminga pulsar (2FGL J0633.9-I-1746) and the Crab (2FGL J0534.5-I-2201). For the W44 count 
map, W44 is the dominant source, sub-dominant, however, to the galactic diffuse emission. 
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Figure 2: (A and B) Gamma-ray spectra of IC 443 (A) and W44 (B) as measured with the 
Fermi-LAI. Color-shaded areas bound by dashed lines denote the best-fit broadband smooth 
broken power law (60 MeV to 2 GeV), gray-shaded bands show systematic errors below 2 
GeV due mainly to imperfect modeling of the galactic diffuse emission. At the high-energy 
end, TeV spectral data points for IC 443 from MAGIC (29) and VERITAS (30) are shown. 
Solid lines denote the best-fit pion-decay gamma-ray spectra, dashed lines denote the best-fit 
bremsstrahlung spectra, and dash-dotted lines denote the best-fit bremsstrahlung spectra when 
including an ad hoc low-energy break at 300 MeV c^^ in the electron spectrum. These fits were 
done to the Fermi LAT data alone (not taking the TeV data points into account). Magenta stars 
denote measurements from the AGILE satellite for these two SNRs, taken from {31) and (79), 
respectively. 
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Figure 3: Proton and gamma-ray spectra determined for IC 443 and W44. Also shown are 
the broadband spectral flux points derived in this study, along with TeV spectral data points for 
IC 443 from MAGIC (29) and VERITAS (30). The curvature evident in the proton distribution 
at ~ 2 GeV is a consequence of the display in energy space (rather than momentum space). 
Gamma-ray spectra from the protons were computed using the energy-dependent cross section 
parameterized by (32). We took into account accelerated nuclei (heavier than protons) as well 
as nuclei in the target gas by applying an enhancement factor of 1.85 (33). Note that models 
of the gamma-ray production via pp interactions have some uncertainty. Relative to the model 
adopted here, an alternative model of (6) predicts ~ 30% less photon flux near 70 MeV; the two 
models agree with each other to better than 15% above 200 MeV. The proton spectra assume 
average gas densities of n = 20 cm~^ (IC 443) and n = 100 cm~^ (W44) and distances of 1.5 
kpc (IC 443) and 2.9 kpc (W44). 
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Analysis of the LAT data 

In its normal mode, the LAT scans the whole sky every 3 hours (two orbits). In this analysis, data 
taken in the sky-survey mode starting from the beginning of scientific operation on 2008 August 
4 to 2012 July 16 were used. These data were analyzed using the LAT Science Tools package 
(v9r29), which is available from the Fermi Science Support Center. The analysis used the P7v6 
version of the instrument response functions which take into account accidental coincidence 
effects in the detectors. Only events passing the "Source" class cuts are used in the analysis, 
and events coming from zenith angles > 100° were discarded to reduce the contribution from 
the emission in the Earth's upper atmosphere. To further reduce these Earth-emitted gamma 
rays, intervals when the Earth was in the field of view were excluded, specifically when the 
rocking angle of the LAT was greater than 52° or when parts of the region-of-interest (ROI) 
were observed at zenith angles > 100°. All gamma rays with energies > 60 MeV within a 
20° X 20° region around the nominal positions of the sources were used. 

During the broad-band fit of the source of interest, all 2FGL sources within the field of view 
were part of the likelihood model with the spectral models as adopted in the 2FGL catalog {14). 
Because the 2FGL catalog corresponds to 2 years of data, while the dataset used in the study 
comprises 4 years, we searched for additional (faint) sources within the ROI by constructing 
maps of residual significance after subtracting all known sources in the field. For the W44 region 
three additional sources were identified at positions (a;j2ooO) '^j2ooo)= (283.92°, 2.79°) with a TS 
value above 60 MeV of 137, (285.25°, 7.11°) with a TS value of 98, and (290.90°, 6.62°) with a 
TS value of 52, where aj2ooo and (5j2ooo are right ascension and declination for the epoch 2000. 
In the IC 443 ROI one additional source was identified at (aj2ooo, 5j2ooo)= (100.17°, 17.80°) with 
a TS value of 1 14. None of these sources affect the spectral parameters fit for the sources under 
study but they were included for completeness and to make the residual maps flat. To account for 
correlations between close-by sources, the normalizations of close-by sources (shown as crosses 
in Figured]) and of the Galactic and isotropic diffuse emission were left free in the broad-band 
and in the bin-by-bin fitting procedure. In the case of IC 443, the Crab was modeled with three 
components for which the normalization was left free following {34): one for the Crab pulsar, 
and two for the synchrotron and for the inverse Compton part of the Crab Nebula respectively. 
The maximum likelihood estimation for spatially and spectrally binned data is performed inside 
a square region of 20° x 20° centered on each SNR (corresponding to the regions shown in 
Figured]). The spatial model for both SNRs was taken as a uniform disk with an angular radius 
of 21' (77). The exact choice of the template has only a small impact on the measurement of 
the spectrum at low energies, due to the large point-spread function of the instrument (the 68% 
containment radius is 6° at 100 MeV), while it affects the spectrum above 1 GeV at the < 30% 
level {11). 



Systematic errors 

For regions near the Galactic equator, the dominant systematic error at energies below ~ 1 GeV 
arises from the uncertainty in the model of the Galactic diffuse emission. To test for this effect 
on the emission assigned to IC 443 and W44, a range of alternative diffuse models calculated by 
GALPROP (35) and then adjusted for agreement with the observed diffuse gamma-ray emission 
within the ROI was used. The GALPROP models were chosen to represent a broad range of 
the parameters scanned in (36). They were adjusted for overall consistency with the all-sky 
gamma-ray data in a likelihood fit. When fitting the alternative diffuse models in this analysis, 
the overall normalization was varied. Eight models were tested: the parameters varied among 
the models are the radial distribution of the cosmic-ray sources (SNR-like or pulsar-like), the 
size of the cosmic -ray halo (4 kpc or 10 kpc) and the spin temperature of atomic hydrogen (150 
K or optically thin). It is not a priori obvious which of these (or the standard) diffuse models 
are the best descriptions of the data for the particular ROIs under study, but the range of diffuse 
models is intended to test systematics related to uncertainties in the model for the Galactic 
diffuse emission. The gray bands in Figure |2] show the systematic error due to uncertainties in 
the Galactic diffuse emission, derived as the envelope of the broad-band fits for the different 
models. 

An additional important source of systematic errors that needs to be considered when ana- 
lyzing LAT data at low energies is the effect of energy dispersion. In the standard analysis of 
LAT data, energy dispersion is typically neglected. As long as the LAT effective collection area 
is relatively independent of energy (above ~ 200 MeV) the bias on the flux in an energy bin that 
arises from misreconstructing the energies of the events is relatively mild (of the order 1-2%). 
However, the combination of energy dispersion and the rapidly changing effective area below 
100 MeV for the instrument response functions used in this study (P7SOURCE_V6) does result 
in biased measurements of flux and spectral index of the source under study. In order to quan- 
tify the effects of energy dispersion, simulations of the IC 443 region were performed taking 
into account the observation conditions of these sources in the real dataset and simulating all 
nearby sources in the region. These simulations include the effect of energy dispersion. The 
analysis performed on these simulations was exactly the same as performed on the real data. 
For the study presented here two effects are important: the low-energy index in the broad-band 
fit (Fi) will be less steeply falling when going to lower energies and the flux points below ~ 150 
MeV will be shifted upwards when ignoring energy dispersion. The simulations can be used to 
quantify the effect when ignoring energy dispersion. For the spectrum of IC 443, Fi is typically 
less steeply falling by ~ 1 (see Figure ISTI) when ignoring the effect. The flux points are higher 
by the following fractions (60 - 81 MeV : 39%, 81 - 110 MeV : 23%, 110 - 150 MeV : 12%). 
For W44 we expect corrections of the same magnitude. At higher energies these factors be- 
come smaller than 5%. These studies demonstrate that taking energy dispersion in the analysis 
into account results in the correct reconstructed values. The magnitude of the bias of the index 
depends on the spectrum of the Galactic diffuse emission model at energies much lower than 
the 60 MeV limit for the present analysis. Unlike the simulations just described we do not have 



a perfect model for the diffuse emission, and so even with energy dispersion included in the 
calculations, Fi is difficult to determine accurately, although for any plausible model the sign 
of the bias in Fi from neglecting energy dispersion is always the same. In order to not overes- 
timate the magnitude of the spectral break we have decided to neglect energy dispersion in the 
spectral analysis for this paper This decision will make the measured low-energy break less 
pronounced but has little impact on the fitted spectrum for the cosmic -ray protons responsible 
for the gamma rays, since it just affects the exact shape of the pion bump. 

Leptonic models 

In addition to 7r°-decay emission, bremsstrahlung and inverse Compton scattering from elec- 
trons are possible gamma-ray emission mechanisms for SNRs. Here we consider the cases in 
which the leptonic processes are dominant in the energy range covered by Ferm«-LAT. 

As already discussed in the main text, it is difficult to explain the large gamma-ray luminos- 
ity of W44 and IC 443 with inverse Compton scattering (see e.g. dashed curve in Figure |S2T). 
Considering typical interstellar radiation fields and infrared photons produced locally by the 
SNR as target photons, the total kinetic energy of electrons is required to be ~ 10^^ erg, which 
means almost 100% of the kinetic energy released by a supernova explosion should be con- 
sumed for electron acceleration. In addition, the measurement of a sharply falling spectrum 
below ~ 200 MeV is inconsistent with the radio spectrum. 

We thus consider models in which electron bremsstrahlung is dominant in the LAT energy 
band. We assume for the electrons a smoothly broken power law, similar to what was used for 
the proton spectrum. The radio data of SNRs IC 443 and W44 strongly constrain the electron 
indices below the break (si). On the other hand, the index above the break (S2) and the location 
of the break pbr can be determined by fitting the gamma-ray spectral shapes. The total number 
of electrons, magnetic field strength, and ambient gas density are then chosen to simultaneously 
fit the synchrotron radio and the gamma-ray bremsstrahlung flux. The magnetic field cannot be 
too low otherwise a break in the synchrotron spectrum corresponding to pbr appears in the radio 
band where pure power-law-type spectra are observed. Figures [S2] and [S3] show the models for 
IC 443 and W44, respectively (dashed lines). An important discrepancy can be found in the 
energy range below < 200 MeV, where we found low-energy breaks in the Fermi-hAI data. 

In order to match the Fermi-hAI data, the electron spectra need to have additional low- 
energy breaks. In Figures |S2] and |S3l we plot bremsstrahlung models with the low-energy 
breaks added to the electron spectra (dash-dotted lines). To make the bremsstrahlung spectra 
below ~ 200 MeV as hard as possible, we have applied abrupt breaks at 300 MeV in 
the electron spectra for both W44 and IC 443. Even with this extreme assumption, the model 
curves do not become as hard as the Fermi-hAI spectra below ~ 200 MeV as seen in Figure |2l 
although the model and data are marginally consistent considering the systematic errors. 

Our results indicate most of the GeV gamma-ray emission is of hadronic origin, i.e., 7r°- 
decay gamma rays (solid lines). The observed synchrotron radio emission still constrains the 
accelerated electrons in the SNRs. As shown in Figures |S2] and [S3l the radio spectra, measured 



at frequencies from several 10 MHz up to 10 GHz, have a power-law shape with spectral indices 
of a = 0.36 and 0.37 for IC 443 and W44, respectively. This means that the radio-emitting elec- 
trons have a power-law index of si = 1.72 and 1 .74 for IC 443 and W44, respectively. The elec- 
tron indices are smaller than the proton indices of si = 2.3-2.4 determined from the gamma-ray 
spectra. It should be noted that the radio-emitting electrons have lower momentum compared 
with the GeV-emitting protons in the crushed cloud scenario of (20, 25), and therefore differ- 
ent indices would be explained by introducing another spectral break. Moreover, momentum 
spectra of electrons and protons do not have the same shape at low energies if re- acceleration 
of the Galactic cosmic rays is responsible for the observed radio and gamma-ray emission (25). 
Given that middle-aged, partially radiative SNRs with molecular cloud interactions are complex 
systems, there is another possibility that the radio-emitting electrons and GeV-emitting protons 
have different origins, particularly in the case of IC 443. Strong radio emission in IC 443 comes 
from the northeastern part of the remnant where the shock is propagating in atomic clouds, 
while the gamma-ray peaks are located near the interacting molecular clouds. 

Energetics 

As discussed in the main text, the dominant sites of gamma-ray emission are likely to be the 
shocked molecular clouds. The mass of the shocked gas is estimated as Mshocked ~ 1 x IO^Mq 
and ~ 5 X IO^Mq for IC 443 and W44, respectively. To explain the gamma-ray luminosity, 
the CR energy density in the clouds should be ~ 400 eV cm""^, much larger than that of the 
Galactic CRs. Therefore, adiabatic compression of the Galactic CRs alone cannot explain the 
required energy density; we need re- acceleration of the pre-existing Galactic CRs (25) and/or 
freshly accelerated particles that have entered the dense region (20). The total energy content 
in protons with p > 0.8 GeV (above the pion production threshold) amounts to Wqr ~ 4 x 
10^*^(71/20 cm-=^)-i erg and 4 x 10^*^(71/100 cm-^)-^ erg for IC 443 and for W44, respectively. 
Here n denotes the effective gas number density, given by n = Mshocked/Kheii where V^heii 
is the volume of the SNR shell that contains shock- accelerated CR particles. We assume that 
the measured gamma-ray emission is produced mainly by CRs with spectra given by Eq. (1) 
interacting with the shocked clouds inside SNRs. Also, the CR density is assumed to be uniform 
in the SNR shell, which contains the shocked clouds. Comparing estimates of the explosion 
energies (Wsn ~ 1 x 10^^ ergs for IC 443 (37) and ~ 5 x 10^^ ergs for W44 (38)) with Wcr, 
we obtain VFcr/W^sn of the order of 1-10% in these objects at their current ages. The ratio 
strongly depends on the (rather uncertain) assumed density at the location of the gamma-ray 
production. 
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Figure S 1 : The effect of ignoring or taking the energy dispersion into account in the fitting of simulated 
data (with energy dispersion) is demonstrated here. For each case 5 realizations of the IC 443 ROI were 
simulated. The results from fitting with the energy dispersion ignored are shown as solid lines, the results 
from fitting with the energy dispersion taken into account are shown as dotted lines. The gray lines show 
the input spectrum. The left panel shows a simulation of IC 443 only, the middle panel a simulation with 
IC 443 and the Galactic and isotropic diffuse and the right panel a simulation with IC 443, Crab and 
Geminga and the Galactic and isotropic diffuse emission. The lower row shows the deviations from the 
input spectium in the reconstructed values of Fi, F2 and the break energy ii^br- The solid points show 
the cases in which the energy dispersion was ignored, the open crosses the cases in which the energy 
dispersion was switched on. As can be seen, the main effect is on Fi as expected which is ~ 1 steeper if 
energy dispersion is taken into account. 
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Figure S2: The spectral energy distribution of IC 443 with a model in which electron bremsstrahlung 
is the dominant radiation process in the gamma-ray band (dashed hne and dash-dotted as in Figure 
2). In the radio band the dominant radiation process is synchrotron emission. The parameters used 
in the calculation for the dashed curve are si = 1.72, S2 = 3.2, = 10 GeV c^^, B = 50 /xG, 
n = 300 cm^ and We = 5x W^"^ erg (> 1 GeV c~^). The dash-dotted curve indicates the same 
model but with an abrupt low-energy break in the electron spectrum at 300 MeVc^^. For comparison 
the best-fit pion-decay model is shown as a solid line. The dotted line shows a combined bremsstrahlung 
and pion-decay model in which Kgp = 0.01 to demonstrate that such a model is consistent with the data. 
Radio data are taken from (39) 
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Figure S3: The same as Figure |S2] but for W44. The parameters are si = 1.74, S2 = 3.7, pbr = 
10 GeV c-^ B = 9Q iiG,n = 650 cm"^, and T^e = 6 x lO^'^ erg (> 1 GeV c^^). Radio data are taken 
from (26). 



